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ABSTRACT 

We report the detection of pulsed gamma-rays from the young, spin-powered radio pulsar PSR 
J2021+3651 using data acquired with the Large Area Telescope (LAT) on the Fermi Gamma-ray 
Space Telescope (formerly GLAST). The light curve consists of two narrow peaks of similar amplitude 
separated by 0.468 ± 0.002 in phase. The first peak lags the maximum of the 2 GHz radio pulse 
by 0.162 ± 0.004 ± 0.01 in phase. The integral gamma-ray photon flux above 100 MeV is (56 ± 
3 ± 11) x 10 _8 cm _2 s _1 . The photon spectrum is well-described by an exponentially cut-off power 
law of the form ^EL = kE^ r e^ E ^ Ec ^ where the energy E is expressed in GeV. The photon index is 
T = 1.5 ± 0.1 ± 0.1 and the exponential cut-off is E c = 2.4 ± 0.3 ± 0.5 GeV. The first uncertainty 
is statistical and the second is systematic. The integral photon flux of the bridge is approximately 
10% of the pulsed emission, and the upper limit on off-pulse gamma-ray emission from a putative 
pulsar wind nebula is < 10% of the pulsed emission at the 95% confidence level. Radio polarization 
measurements yield a rotation measure of RM = 524 ± 4radm -2 but a poorly constrained magnetic 
geometry. Re-analysis of Chandra data enhanced the significance of the weak X-ray pulsations, and 
the first peak is roughly phase-aligned with the first gamma-ray peak. We discuss the emission region 
and beaming geometry based on the shape and spectrum of the gamma-ray light curve combined with 
radio and X-ray measurements, and the implications for the pulsar distance. Gamma-ray emission 
from the polar cap region seems unlikely for this pulsar. 

Subject headings: pulsars: general — pulsars: individual (PSR J2021+3651) — gamma-rays: obser- 
vations 
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1. INTRODUCTION 

The Large Area Telescope (LAT) went into or- 
bit on 2008 June 11 aboard th e Fermi Gamma-ra y 
Space Telescope (formerly GLAST; IAtwood et al. 20091) . 
Gamma-ray pulsatio ns from the pulsar PSR J2021+3651 
(|Roberts et al. 20021 ) were detected during the first 
weeks of commissioning. The LAT was built to ad- 
dress, along with other pressing questions in high en- 
ergy astrophysics, the extent to which gamma-ray emis- 
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sion by pulsars is a rule, and thereby to better under- 
stand the mechanisms by which the kinetic energy of a 
rotating neutron star is transformed into intense beams 
of radiation. The discovery of a pulsed gamma-ray sig- 
nal from PSR J2021+3651 was reported using AGILE, 
the "Astro-rivelatore Gamma a Immagini LEggero" , by 
Halpern et al. (2008). 

Following the detection of at least six high-energy 
gamma-ray pulsars by the Energetic Gamma Ray Ex- 
periment Telescope (EGRET) on the Compton Gamma 
Ray Observatory (see summary by Thompson 2004), 
observers sought positional coincidences of newly- 
discovered pulsars (e.g., Kramer et al. 2003) with 
EGRET catalog sources (Hartman et al. 1999) or 
searched the EGRET error boxes for new radio pulsars 
(e.g., Crawford et al. 2006). 

In one such search, Roberts et al. (2002) discovered 
PSR J2021+3651 with period 103.7 ms within the GeV 
error box associated with 3EG J2021+3716 (Lamb k, Ma- 
comb 1997; Roberts et al. 2001). PSR J2021+3651 is 
young (characteristic timing age 17 kyr) and energetic 
(E = 3.38 x 10 36 ergs s" 1 ). The NE2001 model (Cordes 
& Lazio 2002) for the Galactic distribution of free elec- 
trons for this line of sight (1 = 75?21, b = 0?13) and 
the dispersion measure of DM ps 370 pc cm~ 3 suggest 
a distance of D = 12 kpc, with a fractional uncertainty 
that can exceed 50%. Subsequent investigations with the 
Chandra X-ray Observatory revealed a pulsar wind neb- 
ula (PWN) and possible pulsations in X-rays (Hessels et 
al. 2004). The torus of this "Dragonfly" PWN is clearly 
resolved, providing an estimate of the orientation of the 
pulsar's spin axis relative to the observer's line of sight of 
85° ± 1° (Van Etten et al. 2008). All authors point out 
that such a large distance implies high efficiency r\ for the 
conversion of spin-down power into gamma-rays, becom- 
ing unphysical (r\ > 100%) for some beam scenarios, and 
explore different ways to constrain the distance. The X- 
ray spectra are consistent with a distance of 2 to 4 kpc. 
This paper adds new elements to the discussion: a more 
detailed gamma-ray light curve with phase-resolved spec- 
troscopy, a detailed comparison of the gamma-ray light 
curve with the predictions of various models, and radio 
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polarization measurements. 

PSR J2021+3651 illustrates the importance of sus- 
tained monitoring of pulsar timing. It is one of the more 
than 200 pulsars with E > 10 34 ergs s _1 monitored by 
the program coordinated between the radio and X-ray 
timing community and the Fermi LAT team (Smith et 
al. 2008). PSR J2021+3651, like many young pulsars, 
exhibits timing noise. Not only does contemporaneous 
timing allow accurate comparisons between radio and 
gamma-ray light curves, but folding at a known period 
provides significantly greater sensitivity than searching 
for a periodicity. McLaughlin & Cordes (2004) detected 
significant periodicities for PSR J2021+3651 in two of 
eight EGRET viewing periods (VPs) containing the pul- 
sar by extrapolating the timing solution and searching 
around a range of period and period derivative. The 
light curve was consistent between the two VPs and is 
similar to the LAT light curve. However, they were un- 
able to detect significant pulsations in the other six VPs. 
PSR J2021+3651 was one of only two pulsars detected 
by them in this manner due to the large number of trials 
needed when extrapolating timing solutions back several 
years. The ongoing pulsar monitoring for Fermi greatly 
reduces such problems. 

2. OBSERVATIONS 

The LAT is a pair-production telescope. It consists 
of tungsten foil and silicon microstrip converter/trackers 
(pair conversion and track measurement); hodoscopic ce- 
sium iodide calorimeters (energy measurement); plastic 
scintillator anticoincidence detectors (charged-particle 
rejection); and a programmable trigger and data acquisi- 
tion system. The LAT's excellent sensitivity stems from 
a large effective area (~ 8000 cm 2 ) and superior angular- 
resolution. The broad field of view (2.4 steradian) allows 
long exposures to the whole sky. The LAT is sensitive to 
gamma-rays with energy > 20 MeV. Verification of the 
on-orbit response is continuing but appears consistent 
with expectations 56 . 

Gamma-ray events recorded with the LAT have time- 
stamps that derive from a GPS clock on the Fermi satel- 
lite. Ground tests using cosmic ray muons demonstrated 
that the LAT measures event times with a precision sig- 
nificantly better than 1 /is. On orbit, satellite teleme- 
try indicates comparable accuracy. Transformation to 
the solar system barycenter and phase calculations were 
done with the Fermi LAT "Science Tools", shown to 
be accurate to better than a few /xs for isolated pulsars 
(Smith et al. 2008). Degradation of the barycentered 
time resolution from uncertainty in Fermts position is 
negligible. End-to-end performance of the timing sys- 
tems was confirmed using the bright EGRET pulsars. 
The LAT gamma-ray phases of these pulsars relative to 
the radio phases a gree with previous measurements (e.g., 
iFierro et al. 1998f) . 

Data were acquired in two different Fermi observ- 
ing modes. Between 2008 June 30 and August 3 the 
LAT was often pointed near the northern orbital pole 
{a =18h40m, <5 = 60°, or I = 90°, b = 25°, in mid- July). 
Since then, Fermi has been scanning the entire sky every 
two orbits (approximately 3 hours). The data used here 
were acquired through 2008 November 15. The "diffuse" 



event selection was used ( "Pass 6" , version 1, see Atwood 
et al. 2009), leaving a background due to charged cosmic 
rays comparable to or less than the extragalactic diffuse 
gamma-ray emission. Gamma-rays with measured zenith 
angles greater than 105° were excluded, due to the in- 
tense gamma-ray emission from the Earth's limb caused 
by cosmic rays interacting in the atmosphere. 

PSR J2021+3651 is being observed by the timing con- 
sortium supporting Fermi observations with the NRAO 
Green Bank Telescope (GBT), the NAIC Arecibo tele- 
scope, and the Lovell telescope at Jodrell Bank. The 
combined usage of these observatories provides simul- 
taneously good timing precision and good sampling. 
The most precise timing measurements are acquired at 
Arecibo and GBT, while the best sampling is obtained 
with GBT and Jodrell Bank; the observational setup at 
the latter is described in Hobbs et al. (2004). TheWAPP 
spectrometer used at Arecibo is described in Dowd et al. 
(2000). Here we use a rotational ephemeris based on 
GBT data. 

The phase-connected ephemeris for PSR J2021+3651 
listed in Table 1, contemporaneous with the Fermi obser- 
vations, was derived from 21 observations obtained be- 
tween 2008 June 17 and Nove mber 15 57 . The DE405 so- 
lar system ephemeris was used (jStandish 1 998) . The pul- 
sar was observed at a center frequency of either 1950 MHz 
or 1550 MHz with the GBT Pulsar Spigot (Kaplan et 
al. 2005), yielding total power samples every 81.92 /is 
in each of 768 frequency channels over a bandwidth of 
600 MHz. Each observation lasted for 5 minutes, from 
which we derived a time of arrival with typical uncer- 
tainty of 0.2 ms. PSR J2021+3651 exhibits rotational 
instability that is significant over the span of 5 months. 
We used the TEMPO timing software 58 and describe the 
pulsar rotation well by fitting for its frequency and first 
two derivatives. This timing solution had small unmod- 
eled residual features (xl = 1-4), with a post-fit rms 
of 0.2 ms. The DM is used to correct the time of ar- 
rival of the radio pulse to infinite frequency for absolute 
phase comparison with the gamma-ray profile. Because 
the radio profiles of PSR J2021+3651 are significantly 
scattered at the frequencies that we use for timing, an 
ordinary TEMPO fit biases the DM upwards from its 
true value. In order to correct for this, we fit scatter- 
ing models to the (assumed intrinsic) profile at 5 GHz 
until the resulting profiles matched the observed ones. 
The measured pulse broadening (scaled to 1 GHz) ap- 
pears to be about a factor of two larger than previously 
estimated by Hessels et al. (2004). The corresponding 
DM = 367.5 ± 1 pc cm -3 is 0.5% smaller than the value 
without this correction, where the estimated uncertainty 
is dominated by systematic effects. The resulting phase 
uncertainty after extrapolation from 1.95 GHz to infinite 
frequency is ±0.01. 

In addition to the Spigot observations, we observed 
the pulsar at 2000 MHz using the new Green Bank Ul- 
timate Pulsar Processing Instrument (GUPPI) 59 at the 
GBT for 1.3 hours. GUPPI provides full-polarization 

57 This and other ephemerides used in Fermi results will be avail- 
able from the Fermi Science Support Center (FSSC) data servers 



atlhttp:/ /fermi. gsfc.nasa.gov/ssc 

DB http:/ /www. atnf.csiro.au/research/pulsar/tempo/ 
59 https://wikio.nrao.edu/bin/view/CICADA/GUPPiUsersGuide 
http: //www-glast . slac . Stanford, edu/sof tware/IS/glast_lat_perf ormance .htm 
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spectra in 2048 channels over 800 MHz of bandwidth ev- 
ery 40.96 [is. The resulting data were polarization- and 
flux-calibrated based on observations of a local pulsed 
noise signal and the bright quasar J1445+0958. Analy- 
sis was performed using the PSRCHIVE software package 
(Hotan et al. 2004). 

We re-analysed the 20.8 ks Chandra ACIS-S contin- 
uous clocking data from MJD 52682 described in Hes- 
sels et al. (2004), using the current version of Chan- 
dra analysis software, and a timing solution built from 
near-contemporaneous radio observations acquired with 
the Lovell telescope at Jodrell Bank. The light curve is 
shown in Figure [T] The deviation from a flat light curve 
is slightly more significant (4.5cx) than that reported by 
Hessels et al. (3.7a), has the same overall shape as their 
Figure 3, but is shifted in phase. The first peak seems to 
be within sa 0.1 in phase with the first gamma-ray peak 
albeit with weak statistics. Hessels et al. discuss the 
possibility of a non-thermal component from the X-ray 
point source, and Van Etten et al. measure this compo- 
nent's flux: the ratio of the X-ray non-thermal flux to 
the gamma-ray flux is 10 4 , greater than that observed 
for Vela (10 3 ) or for Geminga (10 2 , Bignami & Caraveo 
1996). The X-ray light curve also resembles those of 
other established gamma-ray pulsars (Kaspi et al. 2006). 

3. RESULTS 

The LAT orientation with respect to the celestial 
sphere has been calibrated to a precision of 30 , using 
more than a dozen known bright point sources. The 
CTiandra-derived position of PSR J2021+3651 is 0162 
from the center of the gamma-ray 68% containment con- 
tour, which has a 0'.75 radius. Figure [2] shows the distri- 
bution of counts for the region surrounding the pulsar. 

3.1. Light curve 

The top frame of Figure Q] shows the phase histogram 
for the gamma-ray events with energies > 100 MeV, 
within an energy-dependent 68% angular containment 
radius of the radio pulsar position. There are two peaks: 
the first peak (PI) is the phase interval 0.13 < 4> < 0.20, 
and P2 is 0.58 < 4> < 0.68. Both peaks are well-fit 
by a Lorentzian function, yielding full widths at half- 
maximum of 0.021 ± 0.002 and 0.053 ± 0.006 rotations, 
respectively, separated by A0 7 = 0.468 ± 0.002. PI lags 
the maximum of the 2 GHz radio peak shown in the bot- 
tom frame by <5</> 7 _ ra dio — 0.162 ± 0.004 ± 0.01 rotations 
of the neutron star. The first uncertainty is statistical, 
the second is due to the uncertainty in the DM. The in- 
termediate frames show three energy bands: 100 MeV to 
1 GeV; 1-3 GeV; and > 3 GeV. PI fades with increasing 
energy, whereas P2 persists. The highest energy photon 
in this sample has 12 GeV, in P2. The peak positions 
are stable with energy to within 0.01 in phase. 

3.2. Spectrum 

The spectral shape, the cutoff energy, and the inte- 
gral energy flux are observables that can be compared 
with pulsar emission models. PSR J2021+3651 is lo- 
cated in the Cygnus region where diffuse emission is 
bright and there are poorly resolved neighboring sources 
(Fig. [2]). While most spectra in the 3rd EGRET catalog 
were well-modeled by a single power-law spectrum, 3EG 



J2021+3716 was a clear exception. A two power-law fit, 
with the break point fixed at 1 GeV, gave a reduced 
X 2 = 0.55, compared to reduced \ 2 = 2-8 for the simple 
power law (Bertsch et al. 2000; Reimer & Bertsch 2001). 
The two spectral indices of 1.23 ± 0.15 and 3.39 ± 0.36 
bracket the 3EG index of 1.86 ± 0.10. The integral pho- 
ton flux of 35 x 10 -8 cm -2 s _1 thus obtained is a bit more 
than half of the 3EG value, presumably because of how 
neighboring sources were han dled. AGILE data c onfirm 
a spectral break near 2 GeV (|Halpern et al. 20081) . 

The different approaches used in analyses of LAT data 
to find the background contribution and to take into ac- 
count the direction- and energy-dependent instrument 
response are outlined in our recent analysis of LAT data 
for Vela (Abdo et al. 2009). One method estimates the 
background and pulsar spectrum by maximizing the joint 
likelihood for the off-pulse data, phases 0.73 < <j) < 1.05, 
and the pulsed data ('on-off'). Another models observed 
neighboring sources and the diffuse Galactic, extragalac- 
tic, and residual charged backgrounds, again in a likeli- 
hood approach ('gtlike', provided with the Fermi science 
tools) . The 'unfolding' method deconvolves the observed 
events from the instrument response. All three methods 
give consistent results. 

Figure [3] shows the result of the likelihood fit to the 
data assuming a power-law spectrum with an exponential 
cutoff, 

— = kE- T e(- E ' E ^ 
dE 

where the energy E is expressed in GeV and with b = 1 , 
using only the on-pulsc data (0.05 < (j) < 0.73). The 
curves are the 'gtlike' results, whereas 'on-off' gives a 
more faithful representation of the individual source data 
points, which are shown. To reduce confusion from the 
diffuse background and neighboring sources, we fit only 
photons above 200 MeV and extrapolate the result to 
100 MeV. The spectral parameters are listed in Table 
[21 including the integral energy flux h — (4.3 ± 0.1) x 
lO-^ergscm^s" 1 for E > 100 MeV. The Pearson x 2 
value for the fit is 1.1. Applying the likelihood ratio 
test against a spectral model with b left free, there is no 
evidence to reject the simple exponential model, a super- 
exponential cutoff (6 = 2) is disfavored at the 3a level, 
and a pure power-law (b = 0) is strongly excluded at 13cr. 

We also measure the spectra of the two pulsar peaks 
individually. Applying the same background model and 
again imposing 6=1 yields the additional curves shown 
in Figure [3] and summarized in Table O P2 appears to 
persist to higher energies in Figure [1] and the fit re- 
sults formally confirm this impression. However, in light 
of the systematic biases at play, the data are also con- 
sistent with no significant spectral difference between 
the two pulses. Bridge emission in the phase interval 
0.26 < 4> < 0.54 is apparent in Figure [TJ It exceeds the 
backgrounds from the diffuse Galactic emission and from 
neighboring sources with ha significance, with a photon 
flux that contributes roughly 10% of the total on-pulse 
flux. Finally, the off-pulse data show no excess above 
the background, allowing us to place an upper limit on 
the flux of a putative gamma-ray PWN of < 10% of the 
phase-averaged emission, at the 95% confidence level. 

The uncertainties in Table 2 and above are statistical. 
Two effects dominate the systematic biases: modeling 
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of the diffuse emission and neighboring sources over the 
several degree radius dictated by the point spread func- 
tion at the lower energy bound ; and uncertainties in 
the energy-dependent effective area. The latter is calcu- 
lated using the test beam-verified Monte Carlo detector 
simulations (Baldini et al. 2007), and verified on-orbit 
using gamma-ray data from the Vela pulsar. For the 
LAT Vela spectral measurements (Abdo et al. 2009), 
the differences between observed and expected on-orbit 
gamma-ray efficiencies led to a uncertainty on the in- 
tegral energy flux of ±Sh/h = 20%. Since then, the 
differences have been found to arise from charge in the 
silicon tracker deposited by cosmic rays in a time window 
around the gamma-ray event. Gamma-ray event recon- 
struction and selection efficiencies below several hundred 
MeV are smaller than predicted, and the fluxes reported 
here will increase for future analyses taking the effect 
into account. The potential bias in the cutoff energy 
is of order ±0.5 GeV and that of the spectral index is 
ST w 0.1. The conclusion that the spectral shape for 
PSR J2021+3651 is most consistent with b = 1 is unaf- 
fected by these issues. 

4. DISCUSSION 

PSR J2021+3651 is among the first gamma-ray pul- 
sars to be studied using Fermi. The gamma-ray observa- 
tions, combined with X-ray images and spectra as well as 
with radio information such as polarization, allow stricter 
comparison with models than was previously possible. 
Improved knowledge of beam geometries will in turn aid 
interpretation of the number counts of radio-loud and 
radio-quiet gamma-ray pulsars that Fermi sees. 

Emission models fall into two classes: Polar cap (PC) 
scenarios (Daugherty & Harding 1996), which place the 
emission very close to the star surface, and outer mag- 
netosphere models. Of the latter, the outer gap (OG) 
picture (Romani & Yadigaroglu 1995; Cheng et al. 2000; 
Hirotani 2005) ascribes the double pulse to emission at 
the boundaries of a single magnetic pole, while the two 
pole caustic model (Dyks & Rudak 2003), a physical re- 
alization of which may be the modern version of the slot- 
gap (SG) picture (Muslimov & Harding 2004; Harding et 
al. 2008), assigns the two pulses to the trailing bound- 
aries of separate magnetic poles. 

Gamma-rays created at the polar caps interact with 
the intense magnetic fields near the neutron star surface, 
resulting in superexponential spectral cut-offs below a 
few GeV, while OG models predict simple exponential 
cutoffs. LAT confirms and refines the spectral break 
measured with EGRET and also seen with AGILE by 
Halpern et al. (2008). The observed emission beyond 10 
GeV and absence of a sharp cut-off indicate outer mag- 
netosphere processes. 

From the phase separation of the peaks in the gamma- 
ray light curve A</> 7 = 0.468, lag from radio A0 7 _ ra( ji O = 
0.16, and X-ray torus-derived viewing angle £ = 85 ± 1°, 
we have bounds on the pulsar geometry that constrain 
the possible emission models. For the PC model A(/> 7 = 
0.5 is natural for large £ and magnetic inclination a; one 
would expect to observe both radio pulses as well, since 
the observer must view the emission at a small angle 
(3 = ( — a to the magnetic pole to intersect the small PC 
beam. Indeed, Figure [Tj shows a faint radio interpulse. 
On the other hand both outer magnetosphere models 



can produce two narrow gamma-ray pulses at the ob- 
served separation for £ = 85°. Narrow gamma-ray pulses 
arise naturally in outer magnetosphere models, where the 
peaks are formed by caustics. For OG and SG, a mag- 
netic inclination of a w 70° is inferred. For this a, the 
second radio pole may be faint ("grazing") or absent. 

An additional radio observation supports the view that 
PSR J2021+3651 is a nearly orthogonal rotator, that is, 
has a large magnetic inclination. Figure [4] shows po- 
larization data from which we determined the rotation 
measure (RM) of J2021+3651 to be 524±4 rad m~ 2 , and 
measured the polarized pulse profile. The data are insuf- 
ficient to const rain the parameters of the Rota ting Vector 
Model (RVM, iRadhakrishnan k. Cooke 19691 ): they are 
compatible with a ~ 70° found above. 

The large lag from the radio pulse is a challenge for sim- 
ple versions of all models using a vacuum dipole magne- 
tosphere. For the PC model we expect A0 7 _ ra( a o = 0, in 
strong disagreement. Both the OG and SC models pre- 
dict A0 7 _ ra rfi o = 0.05 — 0.1, for radio emission at the star 
surface. However, recent work (Johnston & Weisberg 
2006; Karastergiou & Johnston 2007) suggests that ra- 
dio emission from young, Vela-type pulsars occurs over a 
narrow range of altitudes below roughly 100 times the ra- 
dius of the neutron star, whereas for their older brethren, 
it extends to lower altitudes. At high altitudes aberra- 
tion shifts the radio pulse forward and widens the radio 
cone. If only the leading edge of the cone contributes to 
the narrow radio pulse, then the observed lag is easily 
achieved. 

The emission scenarios predict different relations be- 
tween the gamma-ray luminosity L 7 and the spin-down 
power E, making the efficiency 77 = L^/E an additional 
discriminating observable, especially if applied to a large 
sample of gamma-ray pulsars. To obtain L 7 , the pulsar's 
distance D must be known, and the observed integral 
gamma-ray energy flux h must be extrapolated to the 
full sky, that is, some model of the beam shapes must be 
applied. In the past, the lack of geometric constraints 
for the small number of known gamma-ray pulsars led to 
the convention of simply assuming the gamma-ray beam 
swept out a 1 sr solid angle, from which L 7 = hD 2 ; such 
a narrow beam is appropriate to near-surface polar cap 
emission. To better exploit the available data we define, 
following Watters et al. (2008), a correction factor /□ 
along the Earth line-of-sight Ce as 



/ F 7 (a;C,0)sin (C)dCd0 
2/F 7 (a;Cs J 



Jn(a, Ce) = — - r „ , — t — -TJ2 — W 



such that 

L 7 = 4irfahD 2 . (2) 

Fj(a; C, 4>) is the gamma-ray energy flux as a function 
of C and the pulsar rotation phase </>. In the ratio /n, 
the numerator is the total emission over the full sky, and 
the denominator is the phase-averaged flux for the light 
curve seen from Earth. A 1 sr sky coverage corresponds 
to fn = 4^ = 0.08 and isotropic emission gives /n = 1. 
Note that jq > 1 is possible for beams that are narrow in 
4>, extended in £, and/or have average intensity exceeding 
the value sampled at £g. 

Polar cap models tend to have /n < 0.1, while for the 
fan beams of outer magnetosphere models /a is much 
larger. For the OG model we estimate fa ~ 1.05 while 
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the SG model has fa « 1.1 for the observed £ (Watters 
et al. 2008). 

The distance to PSR J2021+3651 is intriguing. The 
NE2001 electron density model assigns a distance of 
12^^° 7 kpc to the large DM of the pulsar, but also greatly 
underestimates the measured scattering timescale, cast- 
ing doubt on this inferred distance. Van Etten et al. 
(2008) place the pulsar at D = 2 to 4 kpc, based on the 
PWN properties and the neutron star thermal emission. 
The large positive RM presented here is consistent with 
a distance at or beyond 4 kpc: of the 9 pulsars within 
~ 10° of this line- of-sight with measured RMs in the 
ATNF database 60 ([Manchester et al. 20051 ) . those with 
D < 5 kpc have negative RMs, while four others with 
D w 6 to 8 kpc have positive RM < 150 rad m -2 and 
DM from 130 to 240 pc cm" 3 . 

For a neutron star moment of inertia of 10 45 g cm we 
obtain 

V = L 7 /E = 0.25 f n (D/4kpc) 2 . (3) 

High-altitude models are preferred by the pulse and spec- 
tral shapes. A 25% efficiency is amongst the largest for 
all known gamma-ray pulsars, taking the /q values as in 
Table 1 of Watters et al. (2008). The nominal DM dis- 
tance imposes small /o, i.e., polar cap emission, whereas 
the smaller distances imposed by the high-altitude mod- 
els leave unexplained the large observed electron column 
density along the line of sight to the Dragonfly nebula. 

The open cluster Berkeley 87 is 0.5° from PSR 
J2021+3651. Prior to the pulsed gamma-ray de- 
tections, it was suspected to be a proton accelera- 
tor that could explain 3EG J2021+3716 and/or 3EG 
J2016+3657, and was s earched for TeV emission (see e.g. 
lAharonian et al. 2006( 1. The LAT localization clearly re- 
futes Berkeley 87 as the dominant gamma-ray emitter in 
this direction, even if the off-pulse emission upper limit is 
near the intensity predicted by Bednarek (2007), leaving 
the door open for future explorations. The detailed maps 
in Schneider et al. (2007) include the position of PSR 
J2021+3651 and their re-examination may allow part of 
the large free electron column density to be accounted 
for. 

5. CONCLUSIONS 

PSR J2021+3651 was detected in gamma-rays in LAT 
data taken during Fermi instrument commissioning, 
with more data accumulated during the first few months 
of the all-sky survey. The extensive radio pulsar timing 
being performed for the Fermi mission facilitated the 
detection and enhanced the quality of the resulting light 
curves as well as their interpretation. Along with the dis- 
covery of a pulsed signal in the gamma-ray source CTA 
1 identifying a previously mysterious source (Abdo et al. 

60 http:/ /www. atnf.csiro.au/research/pulsar/psrcat/ 



2008), this is a good example of the remarkable capabil- 
ity of the LAT to identify Galactic sources. Of the 171 
unidentified sources in the 3rd EGRET catalog, some 
30 are Galactic (|6| < 3°) and have steady fluxes as mea- 
sured by EGRET. The identification of 3EG J2021+3716 
suggests that many of these other sources could also be 
pulsars. 

The high-resolution gamma-ray light curve, the faint 
radio interpulse, and the polarization data, together with 
earlier X-ray images of the torus and jet of the surround- 
ing pulsar wind nebula, allow comparisons with differ- 
ent pulsar models. The rotation measure adds an argu- 
ment in favor of an intermediate pulsar distance. Re- 
analysis of Chandra X-ray data yields an improved light 
curve. Phase-resolved spectral measurements show that 
both peaks cut off exponentially near 2 GeV. Gamma-ray 
emission from the polar cap is the least plausible explana- 
tion at present, even if the outer magnetosphere models 
imply large gamma-ray efficiencies if PSR J2021+3651 is 
indeed more distant than a few kpc. 
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TABLE 1 

Radio ephemeris of PSR J2021+3651 



Pulsar name J2021+3651 

Right ascension (J2000) a 20:21:05.46 

Declination (J2000) a 36:51:04.8 

Pulse frequency (s^ 1 ) 9.6393948581(3) 

Frequency derivative (s~ 2 ) -8.89419(6) X 10~ 1; 

Frequency 2nd derivative (s~ 3 ) 1.09(5) X 10 -21 

Epoch of fit (MJD [TDB]) 54710.0 

Epoch of zero phase reference (MJD [TDB]) b 54715.22791423678 

Range of fit (MJD) 54634.18-54785.9 

Dispersion measure (pc cm~ a ) c 367.5 ± 1 



Note. — The digits in parentheses arc the nominal Icr TEMPO uncertainties. 
a Celestial coordinates are from Chandra observations (Hessels et al. 2004). b This is 
TEMPO'S "TZRMJD" extrapolated to the solar system barycenter at infinite frequency. 
c See § 2 for details of this measurement. 



TABLE 2 

Spectral results for PSR J2021+3651 



Phase interval 




p a 

(10- 8 cm- 2 s- 1 ) 


h b 

(10~ 10 ergs cm~ 2 s _1 ) 


r 


E c c 
(GeV) 


PI 


0.13 to 0.20 


18 ±1 


1.3 ± 0.1 


1.5 ±0.1 


1.9 ±0.3 


P2 


0.58 to 0.68 


21 ± 1 


1.7 ±0.1 


1.5 ±0.1 


2.8 ±0.5 


Bridge 


0.26 to 0.54 


« 10% of Total pulse 








Total pulse 


0.05 to 0.73 


56 ±3 


4.3 ± 0.1 


1.5 ±0.1 


2.4 ±0.3 


Off pulse 


0.73 to 1.05 


< 10% of Total pulse 









a Integral photon flux (> 100 MeV). b Integral energy flux (> 100 McV). Flux systematic biases could lead to 20% increases, see 
tcxt. c Energy of an exponential cut-off to a power-law spectrum with index T. 
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Fig. 1. — Top frame: Light curve for PSR J2021+3651, for gamma-rays with energy > 100 McV within the energy-dependent 68% 
containment radius 6g s % of the pulsar position. Each bin is 0.01 in phase, and 2 rotation cycles are shown. The horizontal dashed line shows 
the average number of counts in the off-pulse phase interval. Three following frames: Light curves in the three indicated energy ranges. 
Second frame from bottom: Phase-aligned Chandra ACIS-S CC X-ray light curve, with the background rate shown by the horizontal 
dashed line. Bottom frame: The upper curve is the 1950 MHz radio profile obtained using the Pulsar Spigot at the Green Bank Telescope. 
The lower curve is the total intensity profile obtained at Arecibo using the WAPP spectrometers with a 3.7 hr integration time and 300 MHz 
of bandwidth centered at 1500 MHz, vertically offset from the GBT curve for clarity. Both curves show evidence for an interpulse at phase 
0.5 with amplitude ~ 5% of the main pulse. 
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Fig. 2. — Gamma-ray counts per square degree, with E > 100 MeV, centered at the Chandra-derived position of PSR J2021+3651. Bin 
sizes vary such that statistical fluctuations are fixed to a signal-to-noise ratio of 10. The cross indicates the C/iandra-derived position of 
PSR J2021+3651. The "X" indicates the position of the open cluster Berkeley 87. The bright object at upper-left is 3EG J2020+4017. 
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Fig. 3.— Spectral energy distribution E 2 dF/dE for PSR J2021+3651 as fit by "gtlikc" assuming a power-law spectrum with an 
exponential cutoff, for PI (dotted), P2 (dot-dash), and Total pulse (solid). The differential values as estimated by "on— off" are given for 
PI (squares), P2 (diamonds), and Total pulse (circles). The error bars are statistical only. For clarity, the points for PI (P2) are plotted 
2.5% lower (higher) in energy than for Total pulse. 
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Fig. 4. — Bottom: Polarization- and flux-calibrated profile for PSR J2021+3651 obtained with the GUPPI spectrometer at the GBT. 
The black trace corresponds to total intensity, the red to linear polarization, and the blue to circular polarization (greatest to least intensity 
at </> S3 0.6, respectively). The phase reference is different from Fig. [T] The faint interpulse is visible near phase fa 0.1. Top: Position 
angle of linear polarization, corrected to the pulsar frame, accounting for the Faraday rotation implied by the large rotation measure of 
RM = 524±4radm" 2 . 



